Open In Colab | Download notebook.

1. Compositacion de Sondajes

Una variable regionalizada \(Z(\mathbf{u})\) puede definirse, no sólo en cada punto del espacio \(\mathbf{u}\), sino que también en una superficie (2D) o en un volumen (3D), \(Z(\mathbf{V})\). La superficie o el volumen sobre el cual se considera la variable regionalizada se denomina soporte. En general, el soporte de las mediciones es muy pequeño (asimilado a un “punto”), mientras que el que interesa en la práctica puede ser más voluminoso (por ejemplo, las unidades selectivas de explotación en evaluación minera o las unidades de remediación en contaminación de suelo). Esta noción es esencial debido a la dependencia que existe entre el soporte y la distribución estadística de los valores, conocida como efecto soporte: los soportes voluminosos presentan una menor cantidad de valores extremos y una mayor cantidad de valores intermedios que los soportes puntuales. Así, la distribución de los valores (en especial, su varianza) depende del soporte sobre el cual está definida la variable regionalizada.

1.1. Instalar PyGeostat

Instalamos pygeostat, un paquete de Python para modelado geoestadístico. Pygeostat está orientado a la preparación de datos espaciales, la automatización de flujos de trabajo geoestadísticos y el modelado utilizando herramientas desarrolladas en el Centre for Computational Geostatistics (CCG)

[ ]:

## La instalación pide re-iniciar el kernel, asi que no asustarse por esto.

!pip install pygeostat

Desafortunadamente, la librería no ha sido actualizada desde el 2015, por lo que debemos descargar y reemplazar en nuestro computador el archivo desurvey.py. Si estas usando Jupyter, deberas encontra este archivo por tu cuenta en tu computador para reemplazarlo a mano (copiar y pegar el nuevo archivo descargado en tu carpeta).

Si estas usando colab, basta con correl el siguiente código:

[ ]:

!wget -q https://raw.githubusercontent.com/arqlm/arqlm.github.io/main/_static/desurvey.py

!mv desurvey.py /usr/local/lib/python3.12/dist-packages/pygeostat/datautils/desurvey.py

1.2. Importar librerias

[1]:

import pygeostat as gs
import numpy as np
import pandas as pd

[2]:

# Leer collares
collar = gs.DataFile(flname='./dataEvYac_collar.csv',readfl=True,x='SURV_X',
                                 y='SURV_Y',z='SURV_Z',dh='HOLEID')
print('Collars:\n',collar.data.head())

# Leer Surveys
survey = gs.DataFile(flname='./dataEvYac_survey.csv',readfl=True,dh='HOLEID')
print('\nSurveys:\n',survey.data.head())

# Leer Assays: puedes usar tu data una vez que hallas encontrado las ug's
rawdata = gs.DataFile(flname='./Data_sin_compositar_con_UGs.csv',readfl=True,dh='dhid',
                               ifrom='from',ito='to')
print('\nRaw data file:\n',rawdata.data.head())

rawdata.data = rawdata.data[rawdata.data['cu_pct'] >= 0] # Filtrar valores negativos

rawdata.data

Collars:
        HOLEID      SURV_X       SURV_Y    SURV_Z    DEPTH
0   CCDDH-001  472560.798  6925645.753  4204.009   800.10
1  CCDDH-002A  472479.288  6925449.980  4289.060  1010.30
2   CCDDH-003  471898.517  6925719.817  4287.854   860.22
3   CCDDH-004  472109.838  6925708.225  4278.548  1100.00
4   CCDDH-005  472233.539  6925602.721  4322.216   918.00

Surveys:
       HOLEID  Along  Azimuth  Inclination
0  CCDDH-001    0.0  209.000      -60.326
1  CCDDH-001   10.0  208.094      -60.150
2  CCDDH-001   20.0  208.210      -60.143
3  CCDDH-001   30.0  208.512      -60.252
4  CCDDH-001   40.0  208.276      -60.215

Raw data file:
        dhid  length  from    to        Este        Norte  Elevación  au_ppm  \
0  98CCD089    12.2   0.0  12.2  472186.686  6925804.447   4220.763  -99.00
1  98CCD089     1.8  12.2  14.0  472187.202  6925805.493   4213.861    0.30
2  98CCD089     2.0  14.0  16.0  472187.343  6925805.770   4211.986    0.47
3  98CCD089     2.0  16.0  18.0  472187.493  6925806.060   4210.013    0.31
4  98CCD089     2.0  18.0  20.0  472187.642  6925806.351   4208.040    0.29

   ag_ppm  cu_pct  aucn_ppm  cucn_ppm Zmin     Alte    Lito   UG
0   -99.0     NaN     -99.0     -99.0  OXI  ILL_CLO  IBX_MM  2.0
1     2.3   0.011     -99.0     -99.0  OXI  ILL_CLO  IBX_MM  2.0
2    16.2   0.032     -99.0     -99.0  OXI  ILL_CLO  IBX_MM  2.0
3     2.3   0.018     -99.0     -99.0  OXI  ILL_CLO  IBX_MM  2.0
4     2.1   0.010     -99.0     -99.0  OXI  ILL_CLO  IBX_MM  2.0

[2]:
dhid length from to Este Norte Elevación au_ppm ag_ppm cu_pct aucn_ppm cucn_ppm Zmin Alte Lito UG
1 98CCD089 1.80 12.2 14.00 472187.202 6925805.493 4213.861 0.30 2.3 0.011 -99.0 -99.0 OXI ILL_CLO IBX_MM 2.0
2 98CCD089 2.00 14.0 16.00 472187.343 6925805.770 4211.986 0.47 16.2 0.032 -99.0 -99.0 OXI ILL_CLO IBX_MM 2.0
3 98CCD089 2.00 16.0 18.00 472187.493 6925806.060 4210.013 0.31 2.3 0.018 -99.0 -99.0 OXI ILL_CLO IBX_MM 2.0
4 98CCD089 2.00 18.0 20.00 472187.642 6925806.351 4208.040 0.29 2.1 0.010 -99.0 -99.0 OXI ILL_CLO IBX_MM 2.0
5 98CCD089 2.00 20.0 22.00 472187.785 6925806.636 4206.066 0.40 2.1 0.010 -99.0 -99.0 OXI ILL_CLO IBX_MM 2.0
... ... ... ... ... ... ... ... ... ... ... ... ... ... ... ... ...
77836 ZVD002 2.00 684.0 686.00 471293.198 6925823.101 3761.674 0.03 1.0 0.002 -99.0 -99.0 BACK PROP ECS 2.0
77837 ZVD002 2.00 686.0 688.00 471293.540 6925824.040 3759.942 0.01 0.7 0.001 -99.0 -99.0 BACK PROP ECS 2.0
77838 ZVD002 2.00 688.0 690.00 471293.882 6925824.980 3758.209 0.01 0.4 0.001 -99.0 -99.0 BACK PROP ECS 2.0
77839 ZVD002 2.00 690.0 692.00 471294.224 6925825.920 3756.477 0.06 1.3 0.004 -99.0 -99.0 BACK PROP ECS 2.0
77840 ZVD002 2.48 692.0 694.48 471294.607 6925826.972 3754.538 0.06 0.5 0.003 -99.0 -99.0 BACK PROP ECS 2.0

76907 rows × 16 columns

Ahora configuramos el largo del compósito final a obtener: en este caso, 6 metros:

[3]:

%%capture
## muteamos la salida para que no sea tan larga (puede borrar %%capture si quieres ver todo)

comps = gs.desurvey.set_comps(rawdata,6.0)  ## Segundo parámetro controla el largo del compósito final

upscaled = gs.desurvey.get_comps(comps,rawdata,{'cu_pct':'continuous','UG':'categorical'}) ## El tercer parámetro indica el tipo de variable (continua o categórica)
rawdata.origdata = rawdata.data.copy() # Guardar una copia de los datos originales por si acaso...
rawdata.data = upscaled # Hacer que los datos ampliados sean el conjunto de datos principal

Ya podemos ver las estadisticas de la data en su nuevo volumen:

[4]:

upscaled_pos = upscaled[upscaled['cu_pct'] >= 0] # Filtrar valores negativos

upscaled_pos.describe()

[4]:
cu_pct UG from to
count 25414.000000 25414.000000 25414.000000 25414.000000
mean 0.145067 2.944519 358.251330 364.251330
std 0.183558 0.999070 296.482491 296.482491
min 0.001000 1.000000 0.000000 6.000000
25% 0.022667 2.000000 108.000000 114.000000
50% 0.086454 3.000000 272.000000 278.000000
75% 0.219667 4.000000 564.000000 570.000000
max 6.723333 4.000000 1470.000000 1476.000000

Solo falta agregar la posicion o centroide a cada muestra compositada:

[5]:

drillholes = gs.desurvey.set_desurvey(collar,survey,'Along','Azimuth','Inclination')


WARNING: drillhole routines are untested with recent Pygeostat changes!
           Proceed with caution!

[6]:

drillhole_items = list(drillholes.items())

[7]:

%%capture
## muteamos la salida para que no sea tan larga (puede borrar %%capture si quieres ver todo)

gs.desurvey.get_desurvey(rawdata,drillholes,x=collar.x,y=collar.y,z=collar.z)
print(rawdata.data.head())

Filtramos algunos algunos tramos de sondaje sin collares:

[8]:

df_filtered = rawdata.data[rawdata.data['SURV_X'] > -100] # Filtrar valores no deseados
df_filtered

[8]:
cu_pct UG dhid from to SURV_X SURV_Y SURV_Z
0 0.020300 2.0 98CCD089 12.2 18.2 472187.343570 6.925806e+06 4211.772172
1 0.009500 2.0 98CCD089 18.2 24.2 472187.740701 6.925807e+06 4205.844155
2 0.006767 2.0 98CCD089 24.2 30.2 472188.089349 6.925807e+06 4199.915106
3 0.010433 2.0 98CCD089 30.2 36.2 472188.470072 6.925808e+06 4193.995578
4 0.011133 2.0 98CCD089 36.2 42.2 472188.878860 6.925809e+06 4188.084383
... ... ... ... ... ... ... ... ...
22343 -999.000000 -999.0 CCDDH-007 558.0 564.0 471959.484334 6.925444e+06 3935.475961
22344 -999.000000 -999.0 CCDDH-007 564.0 570.0 471961.631627 6.925447e+06 3930.686670
22345 -999.000000 -999.0 CCDDH-007 570.0 576.0 471963.789016 6.925450e+06 3925.900536
22346 -999.000000 -999.0 CCDDH-007 576.0 582.0 471965.956499 6.925453e+06 3921.117560
22347 -999.000000 -999.0 CCDDH-007 582.0 588.0 471968.165175 6.925455e+06 3916.340916

16927 rows × 8 columns

Y finalmente, filtramos algunos algunos tramos de sondaje sin ensayos:

[9]:

upscaled_pos = df_filtered[df_filtered['cu_pct'] >= 0] # Filtrar valores negativos

upscaled_pos.describe()

[9]:
cu_pct UG from to SURV_X SURV_Y SURV_Z
count 16882.000000 16882.000000 16882.000000 16882.000000 16882.000000 1.688200e+04 16882.000000
mean 0.165559 3.143822 391.293660 397.293660 472141.581163 6.925472e+06 3966.341566
std 0.196619 0.926312 284.212204 284.212204 302.008871 2.401475e+02 258.580610
min 0.001000 2.000000 0.000000 6.000000 471116.243127 6.924722e+06 2950.021900
25% 0.032333 2.000000 153.190000 159.190000 471927.824777 6.925302e+06 3787.725996
50% 0.116437 4.000000 330.310000 336.310000 472169.814450 6.925478e+06 4011.448263
75% 0.242961 4.000000 597.475000 603.475000 472356.513996 6.925649e+06 4171.006788
max 6.723333 4.000000 1401.050000 1407.050000 473023.885835 6.926165e+06 4445.782649

Es posible ahora visualizar la información compositada en el espacio:

[10]:

import plotly.io as pio
import plotly.graph_objects as go
import plotly.express as px
import plotly.offline as py

[11]:

df_filtered["UG"] = df_filtered["UG"].astype(str)
fig = px.scatter_3d(df_filtered, x='SURV_X', y='SURV_Y', z='SURV_Z', color='UG',category_orders={'UG': sorted(df_filtered.groupby('UG').groups.keys())})
fig.update_traces(marker=dict(size=1.0))

# Mejor estilo para leyenda
fig.update_layout(
    title_text='Unidades Geológicas para Cu',
    legend_title_text='UGs',
    legend=dict(
        itemsizing='constant',
        font=dict(size=14),
        bgcolor='rgba(255,255,255,0.8)',
        bordercolor='black',
        borderwidth=1
    )
)

fig.show()

C:\Users\Alvaro\AppData\Local\Temp\ipykernel_17776\2927694018.py:1: SettingWithCopyWarning:


A value is trying to be set on a copy of a slice from a DataFrame.
Try using .loc[row_indexer,col_indexer] = value instead

See the caveats in the documentation: https://pandas.pydata.org/pandas-docs/stable/user_guide/indexing.html#returning-a-view-versus-a-copy

Data type cannot be displayed: application/vnd.plotly.v1+json

[12]:

continuous = 'cu_pct'

import plotly.io as pio
# pio.renderers.default = "colab"
pio.renderers.default = "notebook" # usar esta linea en jupyter notebook o vscode
import plotly.graph_objects as go
import plotly.express as px
import plotly.offline as py

fig = go.Figure(data=[go.Scatter3d(
    x=df_filtered['SURV_X'],
    y=df_filtered['SURV_Y'],
    z=df_filtered['SURV_Z'],
    marker=dict(color=df_filtered[continuous],
                colorscale=px.colors.sequential.Rainbow[1:],
                cmin=0.0,
                cmax=df_filtered[continuous].quantile(0.95),
                size=1.0, # Esta opcion controla el tamaño de los datos
                colorbar=dict(
                    title='Cu [%]',
                    thickness=20,
                    # Add a border to the colorbar
                    outlinecolor='black',  # Color of the border
                    outlinewidth=2,       # Width of the border
                    bordercolor='white',   # Background border color (if applicable)
                    borderwidth=1         # Background border width
                )
                ),
    mode='markers',
    opacity=1
)])

# otras opciones para controlar el color de los elementos de la figura
fig.update_layout(scene = dict(
                    xaxis = dict(
                        title='Este',
                        backgroundcolor="white",
                        gridcolor="gray",
                        showbackground=True,
                        zerolinecolor="white",),
                    yaxis = dict(
                        title='Norte',
                        backgroundcolor="white",
                        gridcolor="gray",
                        showbackground=True,
                        zerolinecolor="white"),
                    zaxis = dict(
                        title='Elevación',
                        backgroundcolor="white",
                        gridcolor="gray",
                        showbackground=True,
                        zerolinecolor="white",),
                    ),
                    font_family="Times New Roman",
                    font_color="black", hovermode=False
                  )
fig.show()

[13]:

rawdata.data

[13]:
cu_pct UG dhid from to SURV_X SURV_Y SURV_Z
0 0.020300 2.0 98CCD089 12.20 18.20 472187.343570 6.925806e+06 4211.772172
1 0.009500 2.0 98CCD089 18.20 24.20 472187.740701 6.925807e+06 4205.844155
2 0.006767 2.0 98CCD089 24.20 30.20 472188.089349 6.925807e+06 4199.915106
3 0.010433 2.0 98CCD089 30.20 36.20 472188.470072 6.925808e+06 4193.995578
4 0.011133 2.0 98CCD089 36.20 42.20 472188.878860 6.925809e+06 4188.084383
... ... ... ... ... ... ... ... ...
25467 0.006253 2.0 ZVD002 665.28 671.28 -999.000000 -9.990000e+02 -999.000000
25468 0.017613 2.0 ZVD002 671.28 677.28 -999.000000 -9.990000e+02 -999.000000
25469 0.001000 2.0 ZVD002 677.28 683.28 -999.000000 -9.990000e+02 -999.000000
25470 0.001333 2.0 ZVD002 683.28 689.28 -999.000000 -9.990000e+02 -999.000000
25471 0.003108 2.0 ZVD002 689.28 695.28 -999.000000 -9.990000e+02 -999.000000

25472 rows × 8 columns

Ya sabes como exportar la información en un archivo CSV y como hacer un histograma, para verificar los efectos del nuevo soporte en la distribución.

[14]:

upscaled_pos.to_csv('Data_compositada_con_UGs.csv',index=False)